• language pour l’analyse de données
• créé en 1996 par R. Ihaka et R. Gentleman
• développement assuré par des statisticiens rassemblés dans la R Development Core Team
• publié sous licence GNU GPL

Points forts

• multiplateforme
• libre et gratuit
• language interprété
• puissant avec possibilité d’extensions car fonctionne avec un système de packages
• création de graphiques avancés, exportables dans différents formats
• grande communauté d’utilisateurs très actifs

Points faibles

• plus un language de programmation qu’un logiciel, même si des efforts sont faits pour développer de plus en plus d’interfaces graphiques “presse-bouton”
• pas simple d’accès au départ, même si des efforts sont faits pour simplifier les commandes et instructions
• language limité aux statistiques

Installation

• R : https://cran.r-project.org/
• RStudio : https://posit.co/

• Integrated Development Environment (IDE) pour R
• libre, gratuit et multiplateforme

3 zones

• la console (interpréteur R) et le terminal
• les variables d’environnement et l’historique des commandes
• l’aide et la liste des packages

La console

• permet d’interagir avec l’interpréteur R
• les fléches haut et bas du clavier permettent de naviguer dans l’historique des commandes
• touche entrée du clavier pour exécuter la commande tapée

• pour stocker une valeur ou un résultat, il faut créer une variable et lui affecter une valeur
• le signe d’affectation en R est <-

x <- 10

• on affecte la valeur 10 à l’objet (la variable) x
• on appelle x dans la console pour afficher la valeur

x

## [1] 10

• dès que l’on crée une variable elle apparaît dans l’onglet “Environnement”
• les commandes entrées dans la console apparaissent dans l’onglet “History”

• on peut créer autant de variables que l’on souhaite et faire des opérations avec

x <- 10
y <- 100
x + y

## [1] 110

• les noms de variable peuvent contenir des lettres, chiffres et les symboles “.” et “_”
• ne peuvent pas commencer par un chiffre
• ne pas mettre d’accents
• sensible à la casse (x et X sont deux variables différents)

Projet

• lorsque l’on commence un projet, on crée un nouveau projet
• permet de signifier à RStudio le dossier dans lequel tous les fichiers du projet seront stockés (sources et sorties)
• évite les chemins complets/absolus, on peut alors utiliser les chemins relatifs
• pour créer un nouveau projet : File/New Project
• soit créer un nouveau dossier sur l’ordinateur
• soit utiliser un dossier existant

File/New Project

Scripts

• la console c’est pratique mais cela ne conserve pas le code une fois la session fermée
• utilisation de scripts, qui permettent de stocker les commandes, et de les réutiliser
• pour créer un script: File/New File/R Script
• c’est un fichier texte, qu’on enregistre avec l’extension *.R
• pour envoyer le code sélectionné dans la console soit Run soit Ctrl + Entrée
• “#” pour les commentaires

Notebooks

• plus convivial que le script
• mélange de Markdown (html simplifié) et de code (avec sa sortie juste en dessous)
• s’exporte en html ou pdf
• facilite le partage et la lecture du document

File/New File/R Notebook

• on peut changer title qui est le titre du notebook
• ne pas toucher à output
• ce qui est après les trois tirets peut être supprimé, c’est un texte d’accueil standard
• un notebook est composé de texte au format markdown et de blocs de codes, appelés chunks
• insert/R pour insérer un nouveau chunk ou Ctrl + alt + i
• dans le chunk, Ctrl + Entrée pour exécuter la ligne courante et Ctrl + Maj + Entrée pour exécuter tout le chunk

le résultat de l’exécution du chunk est affiché juste en dessous

le rendu en html du notebook, créé automatiquement à chaque enregistrement du fichier

Markdown

• html simplifié
• utilisé pour structurer le Notebook avec des titres, du texte en gras ou italique, des liens, des listes à puces ou numérotées etc.
• RStudio propose un guide rapide : Help/Markdown Quick Reference

Packages

• système d’extension des possibilités offertes par l’installation de base (R Core)
• ces extensions sont développées et maintenues par la communauté
• elles sont regroupées au sein du réseau CRAN (Comprehensive R Archive Netwok)
• pour installer une extension on peut cliquer le bouton Install de l’onglet Packages
• il faut alors taper le nom du package et cliquer sur Install
• ensuite pour pouvoir l’utiliser il faut le charger, avec la commande library(nom_du_package)
• la commande library() est par convention placée dans le premier chunk du notebook (ou en première ligne du script)

Pratique

• créer un nouveau projet
• créer un nouveau notebook
• installer les packages tidyverse et questionr puis les charger

install.packages("tidyverse")
install.packages("questionr")
library(tidyverse)
library(questionr)

Les types de données

• numériques : integer (entier) ou double (décimaux). Par défaut les objets numériques sont considérés comme des double, il faut ajouter un L pour que ce soit un integer
• chaînes de caractères : character
• booléens : logical
• NA : donnée manquante (Not Available)
• typeof(objet) pour connaître le type

a <- "coucou"
print(typeof(a))

## [1] "character"

Le vecteur

• pour stocker plusieurs valeurs dans le même objet (tableau unidimensionnel)
• toutes les valeurs du vecteur doivent être du même type

vec_1 <- c(1, 10, 15)
vec_1

## [1]  1 10 15

c() est une fonction qui concatène les différentes valeurs dans le même vecteur

(vec_2 <- c("coucou", "hello world"))

## [1] "coucou"      "hello world"

Rmq : dans un chunk si on met entre parenthèses la commande ça l’exécute puis affiche la sortie

une opération appliquée à un vecteur s’applique à toutes les valeurs que contient le vecteur

tailles <- c(1.75, 1.54, 1.85, 1.92)
(tailles_cm <- tailles * 100)

## [1] 175 154 185 192

poids <- c(72, 59, 110, 95)
(imc <- poids / tailles ^ 2)

## [1] 23.51020 24.87772 32.14025 25.77040

• possible créer des vecteurs vides, en précisant le type: numeric(), character(), logical(). Si on ajoute un entier en argument cela crée le vecteur avec le nombre d’éléments correspondants, avec des valeurs par défaut (0 pour numeric, “” pour character et FALSE pour logical)

numeric(10)

##  [1] 0 0 0 0 0 0 0 0 0 0

• l’opérateur “:” permet de créer un vecteur d’entiers qui se suivent

1:10

##  [1]  1  2  3  4  5  6  7  8  9 10

pour accéder à un élément du vecteur, utilisation de l’opérateur [] et entre les crochets le numéro de l’élément souhaité (son index dans le vecteur, sachant que la numérotation commence à 1)

tailles

## [1] 1.75 1.54 1.85 1.92

tailles[2]

## [1] 1.54

possible de demander plusieurs éléments via un vecteur

tailles[2:3]

## [1] 1.54 1.85

ou

tailles[c(1, 3)]

## [1] 1.75 1.85

Pratique

A partir des données ci-dessous, calculer le revenu global du foyer puis le revenu par personne.

conjoint1 <- c(1200, 1180, 1750, 2100)
conjoint2 <- c(1450, 1870, 1690, 0)
nb_personnes <- c(4, 2, 3, 2)

Les fonctions

*- pour effectuer des calculs ou accomplir des actions

• une fonction contient un ensemble d’instructions et renvoie un résultat
• on appelle une fonction avec son nom, suivi de parenthèses avec dedans zéro, un ou plusieurs arguments

# la fonction mean calcule la moyenne
mean(tailles)

## [1] 1.765

# la fonction min donne la valeur minimale du vecteur
min(tailles)

## [1] 1.54

Autres fonctions utiles

• length(v): renvoie la longueur de v (nombre d’éléments)
• sum(v): renvoie la somme des valeurs de v
• range(v): renvoie un vecteur avec le min et le max
• unique(v): renvoie un vecteur à partir de v mais en ayant supprimé toutes les valeurs en double
• seq(from, to, by): création d’un vecteur de from à to par intervalles de by
• seq(from, to, length.out): création d’un vecteur de from à to de longueur length.out
• rep(vecteur, nombre de fois)

# création d'une séquence 
seq(1, 20, by=4)

## [1]  1  5  9 13 17

# création d'un vecteur par répétition de vecteur
rep(1:3, 2)

## [1] 1 2 3 1 2 3

(tailles <- c(tailles, NA))

## [1] 1.75 1.54 1.85 1.92   NA

mean(tailles)

## [1] NA

mean(tailles, na.rm=TRUE)

## [1] 1.765

Attention, comme la fonction mean a plusieurs arguments facultatifs, il faut appeler la fonction en utilisant l’argument avec son nom, ici na.rm.

Pratique

A partir des deux vecteurs ci-dessous, calculer la témpérature moyenne et le total des précipitations

temperature <- c(3.4, 4.8, 8.4, 11.4, 15.8, 19.4, 22.2, 21.6, 17.6, 13.4, 7.6, 4.4)
precipitations <- c(47.2, 44.1, 50.4, 74.9, 90.8, 75.6, 63.7, 62, 87.5, 98.6, 81.9, 55.2)

Avec les vecteurs de notes suivants (6 élèves) :

• faire afficher le min et le max dans chaque matière
• calculer la moyenne générale de chaque élève

python <- c(12, 16, 8, 18, 6, 10)
html <- c(14, 9, 13, 15, 17, 11)
bdd <- c(18, 11, 14, 10, 8, 12)

Le data frame

• tableau de données bi-dimensionnel (lignes, colonnes)
• chaque colonne est un vecteur de données
• les colonnes peuvent être de type différent

Prenons comme exemple un jeu de données, hdv2003, fourni par le package questionr chargé précédemment. On charge le jeu de données avec la commande suivante

data("hdv2003")

cela a pour effet de mettre le jeu de données dans l’environnement de travail (chargé en mémoire dans la session courante de R)
hdv2003 est un extrait de l’enquête Histoire de vie réalisée par l’INSEE en 2003. Il contient 2000 individus et 20 variables.

Cliquer sur l’objet dans l’onglet environnement pour voir le tableau de données (on peut aussi appeler la fonction View(hdv2003))

Fonctions utiles

• nouveau nom: df <- hdv2003 - on fait une copie avec comme nom df plutôt que hdv2003 (long à taper)
• nrow(df): renvoie le nombre de lignes
• ncol(df): renvoie le nombre de colonnes
• dim(df): renvoie un vecteur avec les dimensions, cad nrow et ncol
• names(df): renvoie les noms des variables (colonnes)
• str(df): affiche des infos sur le jeu de données - chaque variable avec son type et les premières valeurs (si on clique sur la flèche bleu dans l’environnement cela affiche le contenu de str(df))
• head(df): affiche les premières lignes
• tail(df): affiche les dernières lignes

df <- hdv2003
dim(df)

## [1] 2000   20

names(df)

##  [1] "id"            "age"           "sexe"          "nivetud"      
##  [5] "poids"         "occup"         "qualif"        "freres.soeurs"
##  [9] "clso"          "relig"         "trav.imp"      "trav.satisf"  
## [13] "hard.rock"     "lecture.bd"    "peche.chasse"  "cuisine"      
## [17] "bricol"        "cinema"        "sport"         "heures.tv"

Sélection de données

• pour accéder à une variable : objet$colonne
• renvoie le vecteur correspondant

head(df$age)

## [1] 28 23 59 34 71 35

• on peut aussi se servir du ‘$’ pour créer une nouvelle variable : dans le dataset il y a la variable heures.tv, on va créer la variable minutes.tv

head(df$heures.tv)

## [1] 0 1 0 2 3 2

df$minutes.tv <- df$heures.tv * 60
head(df$minutes.tv)

## [1]   0  60   0 120 180 120

• autrement, l’opérateur [], comme pour les vecteurs, fonctionne, si on veut sélectionner des lignes
• dans cet opérateur il y a maintenant 2 dimensions: les lignes et les colonnes
• si on veut les lignes 10 à 12, et pour ces lignes toutes les colonnes

df[10:12,]

si seulement certaines colonnes

df[10:12, 2:5]

ou par leur nom

df[10:12, c("age", "sexe", "poids")]

Avec critères

• on sélectionne uniquement certaines lignes en fonction d’un ou plusieurs critères
• ci-dessous on veut uniquement les lignes ou la variable age est égale à 25, et pour ces lignes on veut toutes les colonnes

df[df$age == 25, ]

uniquement les données où age == 25 et sexe == femme

df[df$age == 25 & df$sexe == "Femme", ]

• mêmes critères de ligne que le slide précédent, mais on affiche uniquement certaines colonnes

df[df$age == 25 & df$sexe == "Femme", c("age", "sexe", "poids")]

Pratique

A partir du dataset hdv2003

• quel est l’âge moyen ?
• combien y-a-t’il de retraités ?
• combien y-a-t’il de femmes retraitées ?
• quel est l’âge moyen des femmes retraitées ?

Analyse univariée

• c’est l’analyse d’une seule variable
• dépend de son type:
  • quantitatif (comme l’âge ou le revenu par exemple)
  • qualitatif (nombre limité de valeurs possibles, comme le genre ou la profession par exemple)

Variable quantitative

• on va regarder sa distribution avec des indicateurs simples : minimum, maximum, moyenne, médiane, variance, écart-type
• les fonctions sont min(), max(), range(), mean(), median(), var(), sd()
• on peut aussi calculer les quartiles avec la fonction quantile

quantile(df$age)

##   0%  25%  50%  75% 100% 
##   18   35   48   60   97

• on peut directement demander un quartile précis avec l’argument prob

quantile(df$age, prob=0.25)

## 25% 
##  35

• la fonction summary() donne toutes ces informations d’un coup

summary(df$age)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   18.00   35.00   48.00   48.16   60.00   97.00

Pour une meilleure visibilité, un graphique s’impose

• pour visualiser la distribution d’une variable quantitative on peut utiliser un histogramme
• fonction hist(vecteur)

hist(df$age)

• options pour améliorer le graphique : breaks détermine les classes
• on peut passer un nombre (5 par exemple) ou un vecteur (ici [10, 20, 30, … 100])
• include.lowest pour que la valeur inférieure de la première classe soit incluse

hist(df$age, breaks=seq(10, 100, 10), include.lowest=T)

• ajout d’un titre, modification du label des axes et changement de couleur des barres

hist(df$age, breaks=seq(10, 100, 10), include.lowest=T, main="Distribution de l'âge", 
     xlab="Classes d'âge", ylab="Nombre d'observations", col="gray")

Une autre représentation de la distribution d’une variable quantitative est le boxplot (boîte à moustaches), qui reprend les informations des quartiles

• la boîte représente le 1er quartile, la médiane et le 3ème quartile
• donc la boîte représente 50% des observations (il y a 25% des observations en dessous du trait inférieur de la boîte et il y a 25% des observations au dessus du trait supérieur de la boîte)
• les moustaches représentent les valeurs min et max, sachant qu’elles sont au max à 1.5 fois l’inter-quartile. Si une ou plusieurs valeurs sont en dehors des “moustaches” elles sont représentées par des points. On les considère alors souvent comme des outliers (valeurs extrêmes ou aberrantes)

boxplot(df$age)

Variable qualitative

• la fonction table() compte le nombre d’observations de chaque modalité

table(df$sexe)

## 
## Homme Femme 
##   899  1101

table(df$qualif)

## 
##       Ouvrier specialise         Ouvrier qualifie               Technicien 
##                      203                      292                       86 
## Profession intermediaire                    Cadre                  Employe 
##                      160                      260                      594 
##                    Autre 
##                       58

pour comptabiliser les valeurs manquantes il faut ajouter l’argument useNA=“always”

table(df$qualif, useNA="always")

## 
##       Ouvrier specialise         Ouvrier qualifie               Technicien 
##                      203                      292                       86 
## Profession intermediaire                    Cadre                  Employe 
##                      160                      260                      594 
##                    Autre                     <NA> 
##                       58                      347

La fonction summary() appliquée sur une variable qualitative appelle en fait la fonction table()

Si on veut les pourcentages plutôt que les effectifs, on appelle la fonction prop.table() et on lui passe le résultat de l’appelle à la fonction table()

prop.table(table(df$qualif, useNA="always")) * 100

## 
##       Ouvrier specialise         Ouvrier qualifie               Technicien 
##                    10.15                    14.60                     4.30 
## Profession intermediaire                    Cadre                  Employe 
##                     8.00                    13.00                    29.70 
##                    Autre                     <NA> 
##                     2.90                    17.35

ou sinon on procède en 2 étapes

ma_table <- table(df$qualif, useNA="always")
prop.table(ma_table) * 100

## 
##       Ouvrier specialise         Ouvrier qualifie               Technicien 
##                    10.15                    14.60                     4.30 
## Profession intermediaire                    Cadre                  Employe 
##                     8.00                    13.00                    29.70 
##                    Autre                     <NA> 
##                     2.90                    17.35

• si on a chargé la librairie questionr alors il y a la fonction freq()
• la dernière colonne est en ne comptabilisant que les données “valides”, cad sans les NA

freq(df$qualif)

##                            n    % val%
## Ouvrier specialise       203 10.2 12.3
## Ouvrier qualifie         292 14.6 17.7
## Technicien                86  4.3  5.2
## Profession intermediaire 160  8.0  9.7
## Cadre                    260 13.0 15.7
## Employe                  594 29.7 35.9
## Autre                     58  2.9  3.5
## NA                       347 17.3   NA

• freq() à qq paramètres utiles (valid=F pour ne pas afficher la colonne valid, total=T pour afficher une ligne avec le total, et sort=“dec” pour afficher du plus grand au plus petit)

freq(df$qualif, valid=F, total=T, sort="dec")

##                             n     %
## Employe                   594  29.7
## Ouvrier qualifie          292  14.6
## Cadre                     260  13.0
## Ouvrier specialise        203  10.2
## Profession intermediaire  160   8.0
## Technicien                 86   4.3
## Autre                      58   2.9
## NA                        347  17.3
## Total                    2000 100.0

• c’est le graphique en barres (barplot) qui est la représentation graphique de ces effectifs
• la fonction barplot() prend en argument un dataframe issu de la fonction table() (las=3 c’est pour changer l’orientation des labels)

barplot(table(df$qualif), las=3)

• parfois il est préférable de trier les effectifs avant de faire le graphique, avec la fonction sort()
• cex.axis change la taille de police des labels des ticks (cex pour character expansion ratio)

barplot(sort(table(df$qualif)), cex.names=0.5)

Autre façon de représenter les effectifs des différentes modalités, le diagramme de Cleveland, avec la fonction dotchart()

dotchart(sort(table(df$qualif)))

Pratique

• faire une analyse de la variable heures.tv
• faire une analyse de la variable trav.imp
• essayer différents paramètres dans vos graphiques (notamment main, col, xlab, ylab etc.)
• colors() pour afficher la liste de toutes les couleurs de R

Analyse bivariée

• il s’agit d’étudier la relation entre 2 variables
• l’analyse et les graphiques vont dépendre du type des variables

Deux variables qualitatives

• on fait un tableau croisé, avec la fonction table() toujours, mais cette fois en lui passant deux vecteurs

(tb_prof_genre <- table(df$qualif, df$sexe))

##                           
##                            Homme Femme
##   Ouvrier specialise          96   107
##   Ouvrier qualifie           229    63
##   Technicien                  66    20
##   Profession intermediaire    88    72
##   Cadre                      145   115
##   Employe                     96   498
##   Autre                       21    37

• ça a plus de sens en pourcentages, avec la fonction prop.table()
• margin indique l’axe (lignes=1, colonnes=2) par lequel il faut calculer les pourcentages

prop.table(tb_prof_genre, margin=1) * 100

##                           
##                               Homme    Femme
##   Ouvrier specialise       47.29064 52.70936
##   Ouvrier qualifie         78.42466 21.57534
##   Technicien               76.74419 23.25581
##   Profession intermediaire 55.00000 45.00000
##   Cadre                    55.76923 44.23077
##   Employe                  16.16162 83.83838
##   Autre                    36.20690 63.79310

• la fonction round() permet de fixer le nombre de décimales, 2 ci-dessous

round(prop.table(tb_prof_genre, margin=1) * 100, 2)

##                           
##                            Homme Femme
##   Ouvrier specialise       47.29 52.71
##   Ouvrier qualifie         78.42 21.58
##   Technicien               76.74 23.26
##   Profession intermediaire 55.00 45.00
##   Cadre                    55.77 44.23
##   Employe                  16.16 83.84
##   Autre                    36.21 63.79

Changement d’axe pour le calcul: répartition des professions pour chaque genre

round(prop.table(table(df$qualif, df$sexe), margin=2) * 100, 2)

##                           
##                            Homme Femme
##   Ouvrier specialise       12.96 11.73
##   Ouvrier qualifie         30.90  6.91
##   Technicien                8.91  2.19
##   Profession intermediaire 11.88  7.89
##   Cadre                    19.57 12.61
##   Employe                  12.96 54.61
##   Autre                     2.83  4.06

questionr propose 2 fonctions lprop et cprop, avec en plus le total et un formatage par défaut d’une décimale

lprop(tb_prof_genre)

##                           
##                            Homme Femme Total
##   Ouvrier specialise        47.3  52.7 100.0
##   Ouvrier qualifie          78.4  21.6 100.0
##   Technicien                76.7  23.3 100.0
##   Profession intermediaire  55.0  45.0 100.0
##   Cadre                     55.8  44.2 100.0
##   Employe                   16.2  83.8 100.0
##   Autre                     36.2  63.8 100.0
##   Ensemble                  44.8  55.2 100.0

qq paramètres utiles: ajout de %, des effectifs, et changement du nombre de décimales

lprop(tb_prof_genre, percent=T, n=T, digits=2)

##                           
##                            Homme    Femme    Total    n   
##   Ouvrier specialise         47.29%   52.71%  100.00%  203
##   Ouvrier qualifie           78.42%   21.58%  100.00%  292
##   Technicien                 76.74%   23.26%  100.00%   86
##   Profession intermediaire   55.00%   45.00%  100.00%  160
##   Cadre                      55.77%   44.23%  100.00%  260
##   Employe                    16.16%   83.84%  100.00%  594
##   Autre                      36.21%   63.79%  100.00%   58
##   Ensemble                   44.83%   55.17%  100.00% 1653

Test statistique

• ici on peut faire un test d’indépendance du \(\chi^2\) pour savoir si la répartition homme/femme dans chaque métier est proche ou non.
• plus précisément on pose l’hypothèse \(H_0\) que la profession et le genre sont deux variables indépendantes

chisq.test(tb_prof_genre)

## 
##  Pearson's Chi-squared test
## 
## data:  tb_prof_genre
## X-squared = 387.56, df = 6, p-value < 2.2e-16

• X-squared est la valeur de la statistique du \(\chi^2\) pour le tableau donné, c’est une mesure de distance entre le tableau donné et le tableau théorique si les deux variables étaient indépendantes
• df est le nombre de degrés de liberté
• p-value est la probabilité que l’observe \(H_0\). Ici la p-value est inférieure à 5% le seuil standard, donc on rejette \(H_0\), i.e. on en déduit que les deux variables ne sont pas indépendantes.

pour aller plus loin on peut regarder les résidus du test, avec la fonction chisq.residuals() proposée dans questionr (accessible aussi via chisq.test(table(df$qualif, df$sexe))$residuals)

chisq.residuals(tb_prof_genre)

##                           
##                             Homme  Femme
##   Ouvrier specialise         0.52  -0.47
##   Ouvrier qualifie           8.57  -7.73
##   Technicien                 4.42  -3.98
##   Profession intermediaire   1.92  -1.73
##   Cadre                      2.64  -2.38
##   Employe                  -10.43   9.41
##   Autre                     -0.98   0.88

Interprétation

• si la valeur du résidu pour une case est inférieure à -2, alors il y a une sous-représentation de cette case dans le tableau : les effectifs sont significativement plus faibles que ceux attendus sous l’hypothèse d’indépendance
• si le résidu est supérieur à 2, il y a sur-représentatation de cette case
• si le résidu est compris entre -2 et 2, il n’y a pas d’écart à l’indépendance significatif

Représentation graphique du tableau croisé

la fonction mosaicplot()

mosaicplot(tb_prof_genre, las=3, main="Genre / Profession")

On ajoute une couleur en fonction de la valeur du résidu.
Ici en rouge = sous-représentation et en bleu sur-représentation

mosaicplot(tb_prof_genre, shade=T, las=3, main="Genre / Profession")

Pratique

Faire l’analyse de l’importance du travail (variable trav.imp) en fonction du genre

Une variable qualitative et une variable quantitative

• on va chercher à déterminer si les valeurs de la variable quantitative se répartissent différemment ou non selon les modalités de la variable qualitative
• le boxplot permet de voir cela graphiquement: dans l’appel de la fonction on met la variable quantitative puis “~” et ensuite la variable qualitative, donc on représente la dispersion de la variable quantitative en fonction des modalités de la variable qualitative

boxplot(df$age ~ df$sport, ylab="Age", xlab="Pratique d'un sport", col=c("red", "blue"))

on peut ensuite calculer des indicateurs pour chaque modalité

sportifs <- df[df$sport == "Oui",]
cat("Age moyen des sportifs:", mean(sportifs$age))

## Age moyen des sportifs: 40.92531

non_sportifs <- df[df$sport == "Non",]
cat("Age moyen des non-sportifs:", mean(non_sportifs$age))

## Age moyen des non-sportifs: 52.25137

plus simplement avec la fonction tapply()

tapply(df$age, df$sport, mean)

##      Non      Oui 
## 52.25137 40.92531

le premier argument de la fonction est la variable quantitative, le second est la variable qualitative et le troisième la fonction à appliquer. La fonction applique ici “mean()” au vecteur age groupé par modalité.

il est possible de passer une liste de modalites, avec la fonction list() qui prend en argument des clés (étiquettes) et des valeurs (objets)

tapply(df$age, list(genre=df$sexe, sport=df$sport), mean)

##        sport
## genre        Non      Oui
##   Homme 53.33774 40.72629
##   Femme 51.48059 41.13277

les hommes sportifs ont en moyenne 40.73 ans et les sportives 41.13 ans

Test statistique

• on peut tester si l’âge moyen est statistiquement différent selon que l’on soit sportif ou non
• si la variable quantitative suit un loi normale alors on peut faire un test paramétrique de student, sinon il faut faire un test non-paramétrique de Mann-Whitney
• test de normalité : le test de shapiro - \(H_0\) la distribution suit une loi normale

shapiro.test(df$age)

## 
##  Shapiro-Wilk normality test
## 
## data:  df$age
## W = 0.98081, p-value = 9.351e-16

la p-value est inférieure à 5% donc on rejette \(H_0\), il faut utiliser un test non-paramétrique

sinon, on aurait fait le test de Student (\(H_0\): l’âge moyen des deux modalités n’est pas différent)

t.test(df$age ~ df$sport)

## 
##  Welch Two Sample t-test
## 
## data:  df$age by df$sport
## t = 15.503, df = 1600.4, p-value < 2.2e-16
## alternative hypothesis: true difference in means between group Non and group Oui is not equal to 0
## 95 percent confidence interval:
##   9.893117 12.759002
## sample estimates:
## mean in group Non mean in group Oui 
##          52.25137          40.92531

Test de Mann-Whitney : wilcox.test avec argument paired=F car les deux échantillons (sportifs et non-sportifs) ne sont pas les mêmes personnes

wilcox.test(df$age ~ df$sport, paired=F)

## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  df$age by df$sport
## W = 640577, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0

dans les deux tests la p-value est inférieure à 5% donc on rejette \(H_0\), autrement dit l’âge moyen des sportifs est significativement inférieur à celui des non-sportifs

Pratique

Analyser la relation entre le nombre d’heures de télé (heures.tv) et le statut d’occupation (occup)

rmq: dans tapply, après l’appel de la fonction à appliquer il est possible d’ajouter des arguments qui seront transmis à cette fonction à appliquer. Par exemple tapply(x, y, mean, na.rm=T) enverra l’argument na.rm=T à la fonction mean.

Deux variables quantitatives

chargement d’un nouveau jeu de données, proposé par questionr

data("rp2012")

df <- rp2012

Il s’agit de données sur les communues françaises de plus de 2000 habitants

names(df)

##  [1] "code_insee"       "commune"          "code_region"      "region"          
##  [5] "code_departement" "departement"      "pop_tot"          "pop_cl"          
##  [9] "pop_0_14"         "pop_15_29"        "pop_18_24"        "pop_75p"         
## [13] "pop_femmes"       "pop_act_15p"      "pop_chom"         "pop_agric"       
## [17] "pop_indep"        "pop_cadres"       "pop_interm"       "pop_empl"        
## [21] "pop_ouvr"         "pop_scol_18_24"   "pop_non_scol_15p" "pop_dipl_aucun"  
## [25] "pop_dipl_bepc"    "pop_dipl_capbep"  "pop_dipl_bac"     "pop_dipl_bac2"   
## [29] "pop_dipl_sup"     "log_rp"           "log_proprio"      "log_loc"         
## [33] "log_hlm"          "log_sec"          "log_maison"       "log_appart"      
## [37] "age_0_14"         "age_15_29"        "age_75p"          "femmes"          
## [41] "chom"             "agric"            "indep"            "cadres"          
## [45] "interm"           "empl"             "ouvr"             "etud"            
## [49] "dipl_aucun"       "dipl_bepc"        "dipl_capbep"      "dipl_bac"        
## [53] "dipl_bac2"        "dipl_sup"         "resid_sec"        "proprio"         
## [57] "locataire"        "hlm"              "maison"           "appart"

• pour étudier la relation entre deux variables quantitatives, il est pratique de représenter les points sur un graphique
• par exemple, la part de cadres dans la commune et la part de diplômés du supérieur dans la commune, donne

plot(df$cadres, df$dipl_sup, ylab="Diplômés du supérieur", xlab="Cadres")

On observe une relation linéaire entre les deux variables, c’est à dire que lorsque l’une augmente l’autre augmente également, de façon linéraire (pente à peu près fixe)

A l’inverse, si on représente la relation entre la part de cadres et la part de propriétaires, on obtient

plot(df$cadres, df$proprio, xlab="Part de cadres", ylab="Part de propriétaires")

La relation est beaucoup moins évidente

Coefficient de corrélation

Pour déterminer la relation on va calculer le coefficient de corrélation entre les deux variables. Le coefficient de corrélation varie entre -1 et 1, il mesure la force et le sens de la relation entre deux variables. Un coefficient proche de -1 indique une relation négative (lorsque l’une des variables augmente l’autre diminue), et un coefficient proche de 1 indique une relation positive (les deux variables vont dans le même sens).

• si on a une relation linéaire on peut calculer le coefficient de corrélation de Pearson

cor(df$cadres, df$dipl_sup, method="pearson")

## [1] 0.9371629

Le coefficient indique qu’il y a une forte relation linéaire positive.

On peut faire le test statistique correspondant avec la fonction cor.test()
\(H_0\) : le coefficient est égal à 0

cor.test(df$cadres, df$dipl_sup, method="pearson")

## 
##  Pearson's product-moment correlation
## 
## data:  df$cadres and df$dipl_sup
## t = 193.1, df = 5168, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.9337577 0.9403985
## sample estimates:
##       cor 
## 0.9371629

donc on rejette \(H_0\), la corrélation est significativement différente de 0.

Si on calcule la corrélation entre la proportion de cadres et la proportion de propriétaires

cor(df$cadres, df$proprio, method="pearson")

## [1] 0.1622786

elle est nettement plus faible

Si visuellement on constate que la relation n’est pas linéaire (par exemple une relation non-linéaire monotone), alors il est plus opportun que calculer la corrélation selon la méthode des rangs de Spearman

cor(df$cadres, df$dipl_sup, method="spearman")

## [1] 0.9036273

Régression linéaire

Lorsque l’on a une relation linéaire entre deux variables, on peut vouloir estimer la regression linéaire, c’est à dire estimer les paramètres de la droite qui represente au mieux le nuage de points.

lm(df$cadres ~ df$dipl_sup)

## 
## Call:
## lm(formula = df$cadres ~ df$dipl_sup)
## 
## Coefficients:
## (Intercept)  df$dipl_sup  
##      0.9217       1.0816

Donc la constante est égale à 0.9217 et la pente de la droite à 1.0816

Si la sortie de la fonction lm est placée en argument de la fonction summary alors on a davantage d’informations

summary(lm(df$cadres ~ df$dipl_sup))

## 
## Call:
## lm(formula = df$cadres ~ df$dipl_sup)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -33.218  -1.606  -0.172   1.491  13.001 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 0.921661   0.071814   12.83   <2e-16 ***
## df$dipl_sup 1.081636   0.005601  193.10   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.701 on 5168 degrees of freedom
## Multiple R-squared:  0.8783, Adjusted R-squared:  0.8783 
## F-statistic: 3.729e+04 on 1 and 5168 DF,  p-value: < 2.2e-16

pour ajouter cette droite dans le graphique, on utilise la fonction abline()

reg <- lm(df$cadres ~ df$dipl_sup)

plot(df$cadres ~ df$dipl_sup, xlab="Cadres", ylab="Diplômés du supérieur")
abline(reg, col="red")

Pratique

Analyser la relation entre les variables dipl_aucun (proportion de personnes sans diplômes) et proprio

Présentation

• ensemble de librairies qui ont pour objectif de faciliter le chargement, le nettoyage, le tri, l’analyse de données et les représentations graphiques
• il faut installer le package “tidyverse”, ensuite charger la libraire avec la commande library(tidyverse)

ggplot2 (visualisation), dplyr (manipulation des données), tidyr (mise en forme des données), purrr (programmation), readr (importation de données), tibble (tableaux de données), forcats (variables qualitatives), stringr (chaînes de caractères)

Le Tibble

• version moderne du data frame
• les fonctions des librairies du tidyverse acceptent un data frame en entrée, mais retournent un tibble
• par rapport au data frame, le tibble
  • n’a pas de nom de lignes
  • autorise des noms de colonnes avec espaces, caractères spéciaux, nombres etc.
  • a un affichage par défaut sur qq lignes uniquement, avec des informations sur chaque variable
  • affiche un avertissement si on essaie d’accéder à une colonne qui n’existe pas

as_tibble(rp2012)

## # A tibble: 5,170 × 60
##    code_insee commune    code_region region code_departement departement pop_tot
##    <chr>      <chr>      <chr>       <chr>  <chr>            <chr>         <int>
##  1 01004      Ambérieu-… 82          Rhône… 01               Ain           14233
##  2 01007      Ambronay   82          Rhône… 01               Ain            2437
##  3 01014      Arbent     82          Rhône… 01               Ain            3440
##  4 01024      Attignat   82          Rhône… 01               Ain            3110
##  5 01025      Bâgé-la-V… 82          Rhône… 01               Ain            3130
##  6 01027      Balan      82          Rhône… 01               Ain            2785
##  7 01031      Bellignat  82          Rhône… 01               Ain            3642
##  8 01032      Béligneux  82          Rhône… 01               Ain            3172
##  9 01033      Bellegard… 82          Rhône… 01               Ain           11590
## 10 01034      Belley     82          Rhône… 01               Ain            8870
## # ℹ 5,160 more rows
## # ℹ 53 more variables: pop_cl <fct>, pop_0_14 <int>, pop_15_29 <int>,
## #   pop_18_24 <int>, pop_75p <int>, pop_femmes <int>, pop_act_15p <int>,
## #   pop_chom <int>, pop_agric <int>, pop_indep <int>, pop_cadres <int>,
## #   pop_interm <int>, pop_empl <int>, pop_ouvr <int>, pop_scol_18_24 <int>,
## #   pop_non_scol_15p <int>, pop_dipl_aucun <int>, pop_dipl_bepc <int>,
## #   pop_dipl_capbep <int>, pop_dipl_bac <int>, pop_dipl_bac2 <int>, …

un tibble peut également être converti en data frame avec la fonction as.data.frame()

Lecture / écriture d’un fichier de données

• le csv est le format “standard” pour le stockage et le partage de données
• c’est un fichier texte avec un séparateur de colonnes (soit le point-virgule soit la virgule, soit une tabulation)
• ce format est compatible avec tous les logiciels et toutes les plateformes
• il y a deux fonctions pour ouvrir un fichier csv:
  • read_csv : pour les fichiers avec une virgule comme séparateur et le point en décimale
  • read_csv2: pour les fichiers avec le point-virgule comme séparateur et la virgule en décimale (comme en France par exemple)
• si besoin d’options spécifiques, alors on passera par read_delim, qui est la fonction générique
• à noter que RStudio propose une interface graphique pour l’import de données, en cliquant sur le bouton “Import Dataset” dans l’onglet Environment

• pour écrire dans un fichier csv, il y a les fonction write_csv et write_csv2

• parfois on peut vouloir sauvegarder des objets R dans un fichier, pour une utilisation ultérieure ou pour sauver un dataset à un moment donnée. Cela s’appelle la sérialisation, car le fichier va être enregistré au format de R, et il ne sera lisible que par lui. L’extension du fichier sera alors *.Rdata.

• il suffit d’appeler la fonction save, avec comme premier argument le ou les objets à sauvegarder, et ensuite file=“fichier.Rdata”.

• pour charger un fichier Rdata il suffit d’appeler la fonction load(“fichier.Rdata”)

Recodage des variables

• souvent quand on récupère un jeu de données il faut le nettoyer, le recoder, créer de nouvelles variables et le mettre en forme pour l’analyse
• c’est l’activité qui prend le plus de temps lorsque l’on découvre un jeu de données

Renommer les modalités d’une variable qualitative

• on va utiliser la fonction fct_recode (qui vient du package forcats, installé automatiquement lorsque l’on installe tidyverse)

• on va pouvoir renommer et/ou regrouper des modalités

df <- hdv2003
freq(df$qualif)

##                            n    % val%
## Ouvrier specialise       203 10.2 12.3
## Ouvrier qualifie         292 14.6 17.7
## Technicien                86  4.3  5.2
## Profession intermediaire 160  8.0  9.7
## Cadre                    260 13.0 15.7
## Employe                  594 29.7 35.9
## Autre                     58  2.9  3.5
## NA                       347 17.3   NA

Recodage dans le slide suivant

df$qualif_grouped <- fct_recode(
  df$qualif,
  "Ouvrier" = "Ouvrier specialise",
  "Ouvrier" = "Ouvrier qualifie",
  "Interm" = "Technicien",
  "Interm" = "Profession intermediaire"
)
freq(df$qualif_grouped)

##           n    % val%
## Ouvrier 495 24.8 29.9
## Interm  246 12.3 14.9
## Cadre   260 13.0 15.7
## Employe 594 29.7 35.9
## Autre    58  2.9  3.5
## NA      347 17.3   NA

Si on souhaite recoder une modalité en variable manquante il faut lui affecter la valeur NULL

df$qualif_grouped <- fct_recode(df$qualif_grouped, NULL="Autre")
freq(df$qualif_grouped)

##           n    % val%
## Ouvrier 495 24.8 31.0
## Interm  246 12.3 15.4
## Cadre   260 13.0 16.3
## Employe 594 29.7 37.2
## NA      405 20.2   NA

et pour coder les NA en une modalité spécifique, il faut utiliser la fonction fct_explicit_na

df$qualif_grouped <- fct_explicit_na(df$qualif_grouped, na_level="sans_reponse")

## Warning: `fct_explicit_na()` was deprecated in forcats 1.0.0.
## ℹ Please use `fct_na_value_to_level()` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

freq(df$qualif_grouped)

##                n    % val%
## Ouvrier      495 24.8 24.8
## Interm       246 12.3 12.3
## Cadre        260 13.0 13.0
## Employe      594 29.7 29.7
## sans_reponse 405 20.2 20.2

une autre façon pour regrouper les modalités est d’utiliser fct_collapse

df$qualif_rec <- fct_collapse(
  df$qualif, 
  "Ouvrier"=c("Ouvrier specialise", "Ouvrier qualifie"), 
  "Interm"=c("Technicien", "Profession intermediaire"))
freq(df$qualif_rec)

##           n    % val%
## Ouvrier 495 24.8 29.9
## Interm  246 12.3 14.9
## Cadre   260 13.0 15.7
## Employe 594 29.7 35.9
## Autre    58  2.9  3.5
## NA      347 17.3   NA

Rmq :

1. il existe aussi fct_other() pour regroupe un ensemble de modalités dans la catégorie “Other”, et fct_lump() pour regrouper les modalités les moins fréquences dans la catégorie “Other”, cf. documentation de ces deux fonctions.
2. il existe une interface graphique, proposée par questionr, que l’on peut lancer à partir de Addins puis Levels Recoding.

Réordonner les modalités d’une variable qualitative

df$qualif_rec <- fct_relevel(
    hdv2003$qualif,
    "Cadre", "Profession intermediaire", "Technicien", 
    "Employe", "Ouvrier qualifie", "Ouvrier specialise",
    "Autre"
)

freq(df$qualif_rec)

##                            n    % val%
## Cadre                    260 13.0 15.7
## Profession intermediaire 160  8.0  9.7
## Technicien                86  4.3  5.2
## Employe                  594 29.7 35.9
## Ouvrier qualifie         292 14.6 17.7
## Ouvrier specialise       203 10.2 12.3
## Autre                     58  2.9  3.5
## NA                       347 17.3   NA

Egalement possible de réordonner des modalités par rapport à une fonction appliquée à une autre variable
Ci-dessous on réordonne les modalités de la variable \(occup\) en fonction de l’âge moyen.

df$occup_age <- fct_reorder(df$occup, df$age, mean)
boxplot(df$age ~ df$occup_age, xlab="Occupation", ylab="Age", cex.axis=0.65)

Pratique

• dans les données hdv2003, recoder la variable relig afin d’obtenir

                              n    % val%
Pratiquant                  708 35.4 35.4
Appartenance                760 38.0 38.0
Ni croyance ni appartenance 399 20.0 20.0
Rejet                        93  4.7  4.7
NSP                          40  2.0  2.0

• recoder la variable nivetud afin d’obtenir

                                          n    % val%
N'a jamais fait d'etudes                 39  2.0  2.1
Études primaires                        427 21.3 22.6
1er cycle                               204 10.2 10.8
2eme cycle                              183  9.2  9.7
Enseignement technique ou professionnel 594 29.7 31.5
Enseignement superieur                  441 22.0 23.4
NA                                      112  5.6   NA

• ensuite ré-ordonner les modalités de cette variable afin d’obtenir

                                         n    % val%
Enseignement superieur                  441 22.0 23.4
Enseignement technique ou professionnel 594 29.7 31.5
2eme cycle                              183  9.2  9.7
1er cycle                               204 10.2 10.8
Études primaires                        427 21.3 22.6
N'a jamais fait d'etudes                 39  2.0  2.1
NA                                      112  5.6   NA

• trier les modalités de la variable relig selon l’âge median et faire la représentation graphique avec un boxplot

Créer des modalités en fonction des valeurs d’une variable

• l’idée est de créer une nouvelle variable catégorielle en fonction des valeurs d’une ou plusieurs variables existantes
• pour les cas simples on peut utiliser if_else, et pour les cas plus complexes case_when
• if_else prend 3 arguments, le test, la valeur si vrai et la valeur si faux

v <- c(12, 14, 8, 16, 4)
if_else(v > 10, "Supérieur à 10", "Inférieur à 10")

## [1] "Supérieur à 10" "Supérieur à 10" "Inférieur à 10" "Supérieur à 10"
## [5] "Inférieur à 10"

• dans le test on peut mettre plusieurs conditions

df$statut <- if_else(
    df$sexe == "Homme" & df$age > 60,
    "Homme de plus de 60 ans",
    "Autre"
)

freq(df$statut)

##                            n    % val%
## Autre                   1778 88.9 88.9
## Homme de plus de 60 ans  222 11.1 11.1

case_when étend if_else à plusieurs conditions, sous la forme condition ~ valeur, les tests sont effectués dans l’ordre, donc il faut aller du plus spécifique au plus large.

df$statut <- case_when(
    df$age > 60 & df$sexe == "Homme" ~ "Homme de plus de 60 ans",
    df$age > 60 & df$sexe == "Femme" ~ "Femme de plus de 60 ans",
    TRUE ~ "Autre"
)

freq(df$statut)

##                            n    % val%
## Autre                   1512 75.6 75.6
## Femme de plus de 60 ans  266 13.3 13.3
## Homme de plus de 60 ans  222 11.1 11.1

Comme TRUE est toujours vrai toutes les autres combinaisons sont codées “Autre”.

Découper une variable numérique en classes

• par exemple transformer une variable qui représente le revenu mensuel en une variable qualitative avec des catégories (moins de 500€, de 500 à 1000€ etc.)
• utilisation de la fonction cut
• l’argument breaks va permettre de définir les classes

soit le nombre avec breaks=nombre

hdv2003$age_classe <- cut(hdv2003$age, breaks=5)
freq(hdv2003$age_classe)

##               n    % val%
## (17.9,33.8] 454 22.7 22.7
## (33.8,49.6] 628 31.4 31.4
## (49.6,65.4] 556 27.8 27.8
## (65.4,81.2] 319 16.0 16.0
## (81.2,97.1]  43  2.1  2.1

soit les bornes avec un vecteur

hdv2003$age_classe2 <- cut(
    hdv2003$age, 
    breaks = c(18, 25, 35, 45, 55, 65, 100), 
    include.lowest = T
)
freq(hdv2003$age_classe2)

##            n    % val%
## [18,25]  191  9.6  9.6
## (25,35]  338 16.9 16.9
## (35,45]  390 19.5 19.5
## (45,55]  414 20.7 20.7
## (55,65]  305 15.2 15.2
## (65,100] 362 18.1 18.1

Rmq: le package questionr a développé une interface graphique pour la fonction cut: Addins/Numeric range dividing

Pratique

• avec if_else, créer la variable cinema_bd qui groupe les personnes qui vont au cinéma et qui lisent des bd. Les autres sont mis dans “Autre”. Le résultat attendu est

               n    % val%
Autre        1971 98.6 98.6
Cinéma et BD   29  1.5  1.5

• avec case_when créer la variable genre_fratrie ayant les modalités : Homme avec plus de 2 frères et soeurs, Femme avec plus de 2 frères et soeurs, Autre. Le résulat attendu est

                                         n    % val%
Autre                                  1001 50.0 50.0
Femme avec plus de 2 frères et soeurs  546  27.3 27.3
Homme avec plus de 2 frères et soeurs  453  22.7 22.7

• avec case_when, créer une nouvelle variable avec les modalités : Homme de plus de 30 ans, Homme de plus de 40 ans satisfait par son travail, Femme pratiquant le sport ou le bricolage, Autre. Le résultat attendu est

                                                   n    % val%
Autre                                             714 35.7 35.7
Femme pratiquant le sport ou le bricolage         549 27.5 27.5
Homme de plus de 30 ans                           610 30.5 30.5
Homme de plus de 40 ans satisfait par son travail 127  6.3  6.3

• découper la variable heures.tv en classes afin d’obtenir le tableau de fréquences suivant

       n   %    val%
[0,1]  684 34.2 34.3
(1,2]  535 26.8 26.8
(2,4]  594 29.7 29.8
(4,6]  138  6.9  6.9
(6,12]  44  2.2  2.2
NA       5  0.2   NA

Mise en forme des données

pour pouvoir manipuler les données il faut qu’elles soient bien rangées, ie “tidy”, notamment:

• chaque variable est dans une colonne unique
• chaque colonne contient une unique variable
• chaque ligne correspond à une observation pour chaque variable
• les cellules du tableau représentent les valeurs de chaque observation pour chaque variable

exemple de données non tidy : Population de 3 pays sur 4 années

pour que ce soit tidy il faudrait

• une observation par ligne (ici 4)
• une variable année et une variable population

Le dataset serait alors

Fonction pivot_longer pour créer un dataset “plus long” c’est à dire avec plus de lignes

##   country     2002     2007
## 1 Belgium 10311970 10392226
## 2  France 59925035 61083916
## 3 Germany 82350671 82400996

les colonnes 2002 et 2007 devraient être dans une seule colonne “année” et les valeurs dans une colonne “population”

pivot_longer(data=tidy_ex, cols=c("2002", "2007"), names_to="annee", values_to="population")

## # A tibble: 6 × 3
##   country annee population
##   <chr>   <chr>      <dbl>
## 1 Belgium 2002    10311970
## 2 Belgium 2007    10392226
## 3 France  2002    59925035
## 4 France  2007    61083916
## 5 Germany 2002    82350671
## 6 Germany 2007    82400996

Fonction pivot_wider pour à l’inverse créer un dataset “plus large”, c’est à dire avec plus de colonnes

##   country continent year variable        value
## 1 Belgium    Europe 2002  lifeExp       78.320
## 2 Belgium    Europe 2002      pop 10311970.000
## 3 Belgium    Europe 2007  lifeExp       79.441
## 4 Belgium    Europe 2007      pop 10392226.000
## 5  France    Europe 2002  lifeExp       79.590
## 6  France    Europe 2002      pop 59925035.000
## 7  France    Europe 2007  lifeExp       80.657
## 8  France    Europe 2007      pop 61083916.000

la colonne variable devrait être divisée en deux variables lifeExp et pop

pivot_wider(data=tidy_ex1, names_from=variable, values_from=value)

## # A tibble: 4 × 5
##   country continent  year lifeExp      pop
##   <chr>   <chr>     <dbl>   <dbl>    <dbl>
## 1 Belgium Europe     2002    78.3 10311970
## 2 Belgium Europe     2007    79.4 10392226
## 3 France  Europe     2002    79.6 59925035
## 4 France  Europe     2007    80.7 61083916

Fonction separate pour séparer une colonne en plusieurs

##               eleve  note
## 1 Christophe Martin 10/20
## 2    Jacques Dupont 12/20
## 3 Stéphane Nimporte  8/20

Séparation de élève en prenom et nom

(tidy_ex2 <- separate(tidy_ex2, col=eleve, into=c("prenom", "nom"), sep=" "))

##       prenom      nom  note
## 1 Christophe   Martin 10/20
## 2    Jacques   Dupont 12/20
## 3   Stéphane Nimporte  8/20

Séparation de note en note et note_sur

(tidy_ex2 <- separate(tidy_ex2, col="note", into=c("note", "note_sur"), sep="/"))

##       prenom      nom note note_sur
## 1 Christophe   Martin   10       20
## 2    Jacques   Dupont   12       20
## 3   Stéphane Nimporte    8       20

Fonction separate_rows pour séparer une colonne en plusieurs lignes

##       prenom      nom classe   notes
## 1 Christophe   Martin   6ème 5,16,11
## 2    Jacques   Dupont   5ème      15
## 3   Stéphane Nimporte   4ème   11,17

(tidy_ex3 <- separate_rows(tidy_ex3, notes, sep=","))

## # A tibble: 6 × 4
##   prenom     nom      classe notes
##   <chr>      <chr>    <chr>  <chr>
## 1 Christophe Martin   6ème   5    
## 2 Christophe Martin   6ème   16   
## 3 Christophe Martin   6ème   11   
## 4 Jacques    Dupont   5ème   15   
## 5 Stéphane   Nimporte 4ème   11   
## 6 Stéphane   Nimporte 4ème   17

Fonction unite pour faire l’inverse de separate

tidy_ex3_bis <- tidy_ex3
unite(tidy_ex3_bis, col=eleve, prenom, nom, sep=" ")

## # A tibble: 6 × 3
##   eleve             classe notes
##   <chr>             <chr>  <chr>
## 1 Christophe Martin 6ème   5    
## 2 Christophe Martin 6ème   16   
## 3 Christophe Martin 6ème   11   
## 4 Jacques Dupont    5ème   15   
## 5 Stéphane Nimporte 4ème   11   
## 6 Stéphane Nimporte 4ème   17

Fonction complete pour compléter les combinaisons manquantes

## # A tibble: 6 × 4
##   prenom     nom      matiere note 
##   <chr>      <chr>    <chr>   <chr>
## 1 Christophe Martin   Python  5    
## 2 Christophe Martin   R       16   
## 3 Christophe Martin   HTML    11   
## 4 Jacques    Dupont   HTML    15   
## 5 Stéphane   Nimporte R       11   
## 6 Stéphane   Nimporte Python  17

Jacques n’a pas de note en R et Python, et Stéphane n’a pas de note en HTML

(tidy_ex4 <- complete(tidy_ex4, nesting(prenom, nom), matiere))

## # A tibble: 9 × 4
##   prenom     nom      matiere note 
##   <chr>      <chr>    <chr>   <chr>
## 1 Christophe Martin   HTML    11   
## 2 Christophe Martin   Python  5    
## 3 Christophe Martin   R       16   
## 4 Jacques    Dupont   HTML    15   
## 5 Jacques    Dupont   Python  <NA> 
## 6 Jacques    Dupont   R       <NA> 
## 7 Stéphane   Nimporte HTML    <NA> 
## 8 Stéphane   Nimporte Python  17   
## 9 Stéphane   Nimporte R       11

l’argument nesting permet de signifier que prenom et nom identifient une même personne, sinon il y aurait toutes les combinaisons de prenom et nom en plus dans le dataset

Manipuler les données

le package dplyr propose 5 verbes qui facilitent la manipulation des données (tri, sélection etc.) : slice, filter, select, arrange et mutate

Pour les exemples, utilisation des datasets du package nycflights13 (données des vols de départ de 3 aéroports de New-York en 2013)
il faut installer le package “nycflights13” et le charger

data("flights")
data("airports")
data("airlines")

slice

Permet de sélectionner des lignes du dataset par leur index

slice(airports, 2:8)

## # A tibble: 7 × 8
##   faa   name                                  lat    lon   alt    tz dst   tzone
##   <chr> <chr>                               <dbl>  <dbl> <dbl> <dbl> <chr> <chr>
## 1 AAP   Andrau Airpark                       29.7  -95.6    79    -6 A     Amer…
## 2 ABE   Lehigh Valley International Airport  40.7  -75.4   393    -5 A     Amer…
## 3 ABI   Abilene Regional Airport             32.4  -99.7  1791    -6 A     Amer…
## 4 ABL   Ambler Airport                       67.1 -158.    334    -9 A     Amer…
## 5 ABQ   Albuquerque International Sunport    35.0 -107.   5355    -7 A     Amer…
## 6 ABR   Aberdeen Regional Airport            45.4  -98.4  1302    -6 A     Amer…
## 7 ABY   Southwest Georgia Regional Airport   31.5  -84.2   197    -5 A     Amer…

filter

Permet de sélectionner des lignes du dataset en fonction de critères
Fonctionne avec du boolean indexing, c’est à dire renvoie les lignes pour lesquelles là (ou les) conditions sont vraies

Exemple: tous les vols du mois de janvier

filter(flights, month==1)

## # A tibble: 36,020 × 19
##     year month   day dep_time sched_dep_time dep_delay arr_time sched_arr_time
##    <int> <int> <int>    <int>          <int>     <dbl>    <int>          <int>
##  1  2023     1     1        1           2038       203      328              3
##  2  2023     1     1       18           2300        78      228            135
##  3  2023     1     1       31           2344        47      500            426
##  4  2023     1     1       33           2140       173      238           2352
##  5  2023     1     1       36           2048       228      223           2252
##  6  2023     1     1      503            500         3      808            815
##  7  2023     1     1      520            510        10      948            949
##  8  2023     1     1      524            530        -6      645            710
##  9  2023     1     1      537            520        17      926            818
## 10  2023     1     1      547            545         2      845            852
## # ℹ 36,010 more rows
## # ℹ 11 more variables: arr_delay <dbl>, carrier <chr>, flight <int>,
## #   tailnum <chr>, origin <chr>, dest <chr>, air_time <dbl>, distance <dbl>,
## #   hour <dbl>, minute <dbl>, time_hour <dttm>

Plusieurs conditions

Exemple: tous les vols qui ont entre 10 et 15 minutes de retard au décollage

filter(flights, dep_delay >= 10 & dep_delay <= 15)

## # A tibble: 17,829 × 19
##     year month   day dep_time sched_dep_time dep_delay arr_time sched_arr_time
##    <int> <int> <int>    <int>          <int>     <dbl>    <int>          <int>
##  1  2023     1     1      520            510        10      948            949
##  2  2023     1     1      613            600        13      834            836
##  3  2023     1     1      810            800        10     1153           1147
##  4  2023     1     1      811            800        11     1122           1136
##  5  2023     1     1      814            800        14     1048           1058
##  6  2023     1     1      833            820        13     1436           1435
##  7  2023     1     1      840            828        12     1136           1140
##  8  2023     1     1      858            845        13     1042           1045
##  9  2023     1     1      859            844        15     1207           1209
## 10  2023     1     1     1005            950        15     1226           1235
## # ℹ 17,819 more rows
## # ℹ 11 more variables: arr_delay <dbl>, carrier <chr>, flight <int>,
## #   tailnum <chr>, origin <chr>, dest <chr>, air_time <dbl>, distance <dbl>,
## #   hour <dbl>, minute <dbl>, time_hour <dttm>

Possible de passer le résultat d’une fonction dans la condition

Exemple: les vols avec la distance la plus longue

filter(flights, distance == max(distance))

## # A tibble: 667 × 19
##     year month   day dep_time sched_dep_time dep_delay arr_time sched_arr_time
##    <int> <int> <int>    <int>          <int>     <dbl>    <int>          <int>
##  1  2023     1     1      949            900        49       NA           1525
##  2  2023     1     1     1023           1000        23     1637           1610
##  3  2023     1     2      919            900        19     1543           1525
##  4  2023     1     2      951           1000        -9     1620           1610
##  5  2023     1     3      922            900        22     1535           1525
##  6  2023     1     3     1007           1000         7     1630           1610
##  7  2023     1     4      912            900        12     1511           1525
##  8  2023     1     4     1001           1000         1     1630           1610
##  9  2023     1     5      854            900        -6     1454           1525
## 10  2023     1     5      949           1000       -11     1600           1610
## # ℹ 657 more rows
## # ℹ 11 more variables: arr_delay <dbl>, carrier <chr>, flight <int>,
## #   tailnum <chr>, origin <chr>, dest <chr>, air_time <dbl>, distance <dbl>,
## #   hour <dbl>, minute <dbl>, time_hour <dttm>

select

Permet de sélectionner les colonnes du dataset

Exemple: sélectionner des colonnes latitude et longitude

select(airports, lat, lon)

## # A tibble: 1,251 × 2
##      lat    lon
##    <dbl>  <dbl>
##  1  29.7  -85.0
##  2  29.7  -95.6
##  3  40.7  -75.4
##  4  32.4  -99.7
##  5  67.1 -158. 
##  6  35.0 -107. 
##  7  45.4  -98.4
##  8  31.5  -84.2
##  9  41.3  -70.1
## 10  31.6  -97.2
## # ℹ 1,241 more rows

Si on fait précéder le nom de la colonne d’un “-” cela affiche toutes les colonnes sauf elle

select(airports, -lat)

## # A tibble: 1,251 × 7
##    faa   name                                   lon   alt    tz dst   tzone     
##    <chr> <chr>                                <dbl> <dbl> <dbl> <chr> <chr>     
##  1 AAF   Apalachicola Regional Airport        -85.0    20    -5 A     America/N…
##  2 AAP   Andrau Airpark                       -95.6    79    -6 A     America/C…
##  3 ABE   Lehigh Valley International Airport  -75.4   393    -5 A     America/N…
##  4 ABI   Abilene Regional Airport             -99.7  1791    -6 A     America/C…
##  5 ABL   Ambler Airport                      -158.    334    -9 A     America/A…
##  6 ABQ   Albuquerque International Sunport   -107.   5355    -7 A     America/D…
##  7 ABR   Aberdeen Regional Airport            -98.4  1302    -6 A     America/C…
##  8 ABY   Southwest Georgia Regional Airport   -84.2   197    -5 A     America/N…
##  9 ACK   Nantucket Memorial Airport           -70.1    47    -5 A     America/N…
## 10 ACT   Waco Regional Airport                -97.2   516    -6 A     America/C…
## # ℹ 1,241 more rows

La fonction select peut prendre comme argument des éléments permettant une sélection selon des critères

select(flights, starts_with("dep_"))

## # A tibble: 435,352 × 2
##    dep_time dep_delay
##       <int>     <dbl>
##  1        1       203
##  2       18        78
##  3       31        47
##  4       33       173
##  5       36       228
##  6      503         3
##  7      520        10
##  8      524        -6
##  9      537        17
## 10      547         2
## # ℹ 435,342 more rows

il y a aussi ends_with, contains et matches

Pour renommer les colonnes il y a la fonction rename

rename(airports, longitude=lon, latitude=lat, altitude=alt)

## # A tibble: 1,251 × 8
##    faa   name                      latitude longitude altitude    tz dst   tzone
##    <chr> <chr>                        <dbl>     <dbl>    <dbl> <dbl> <chr> <chr>
##  1 AAF   Apalachicola Regional Ai…     29.7     -85.0       20    -5 A     Amer…
##  2 AAP   Andrau Airpark                29.7     -95.6       79    -6 A     Amer…
##  3 ABE   Lehigh Valley Internatio…     40.7     -75.4      393    -5 A     Amer…
##  4 ABI   Abilene Regional Airport      32.4     -99.7     1791    -6 A     Amer…
##  5 ABL   Ambler Airport                67.1    -158.       334    -9 A     Amer…
##  6 ABQ   Albuquerque Internationa…     35.0    -107.      5355    -7 A     Amer…
##  7 ABR   Aberdeen Regional Airport     45.4     -98.4     1302    -6 A     Amer…
##  8 ABY   Southwest Georgia Region…     31.5     -84.2      197    -5 A     Amer…
##  9 ACK   Nantucket Memorial Airpo…     41.3     -70.1       47    -5 A     Amer…
## 10 ACT   Waco Regional Airport         31.6     -97.2      516    -6 A     Amer…
## # ℹ 1,241 more rows

nouveau nom = ancien nom
si le nouveau nom (ou l’ancien) comporte un ou des espace(s) il faut mettre entre guillemets

arrange

Pour ordonner le dataset selon les valeurs d’une ou plusieurs colonnes

arrange(flights, dep_delay)

## # A tibble: 435,352 × 19
##     year month   day dep_time sched_dep_time dep_delay arr_time sched_arr_time
##    <int> <int> <int>    <int>          <int>     <dbl>    <int>          <int>
##  1  2023     8    16     1839           1929       -50     2254           2218
##  2  2023    12     2     2147           2225       -38     2332             35
##  3  2023     9     3     1801           1834       -33     1934           2040
##  4  2023     8    13     1804           1834       -30     1947           2040
##  5  2023     1    28     1930           1959       -29     2308           2336
##  6  2023     1    18     2103           2129       -26     2301           2336
##  7  2023     5     5     1129           1155       -26     1337           1420
##  8  2023     4    16      925            950       -25     1102           1140
##  9  2023    11    14     2055           2120       -25     2215           2239
## 10  2023    12     6     2234           2259       -25     2350             30
## # ℹ 435,342 more rows
## # ℹ 11 more variables: arr_delay <dbl>, carrier <chr>, flight <int>,
## #   tailnum <chr>, origin <chr>, dest <chr>, air_time <dbl>, distance <dbl>,
## #   hour <dbl>, minute <dbl>, time_hour <dttm>

Utilisation de la fonction desc() pour inverser l’ordre (décroissant)

arrange(flights, desc(dep_delay))

## # A tibble: 435,352 × 19
##     year month   day dep_time sched_dep_time dep_delay arr_time sched_arr_time
##    <int> <int> <int>    <int>          <int>     <dbl>    <int>          <int>
##  1  2023    12    17     1953           1340      1813     2155           1543
##  2  2023    10     1     1240            659      1781     1407            835
##  3  2023     4    25     1201            659      1742     1315            818
##  4  2023     2     7     2045           1700      1665     2352           2025
##  5  2023     4    20      926            619      1627     1135            822
##  6  2023    10    29      856            600      1616     1050            805
##  7  2023     4    30     1818           1617      1561     2001           1820
##  8  2023     3    17     2027           1830      1557     2346           2139
##  9  2023     3    29      633            525      1508      820            700
## 10  2023     3    19     1154           1201      1433     1321           1335
## # ℹ 435,342 more rows
## # ℹ 11 more variables: arr_delay <dbl>, carrier <chr>, flight <int>,
## #   tailnum <chr>, origin <chr>, dest <chr>, air_time <dbl>, distance <dbl>,
## #   hour <dbl>, minute <dbl>, time_hour <dttm>

Plusieurs critères

arrange(flights, month, (dep_delay))

## # A tibble: 435,352 × 19
##     year month   day dep_time sched_dep_time dep_delay arr_time sched_arr_time
##    <int> <int> <int>    <int>          <int>     <dbl>    <int>          <int>
##  1  2023     1    28     1930           1959       -29     2308           2336
##  2  2023     1    18     2103           2129       -26     2301           2336
##  3  2023     1    30     2105           2129       -24     2340           2350
##  4  2023     1    25     1844           1905       -21     2200           2222
##  5  2023     1    27      539            600       -21      700            720
##  6  2023     1    19      820            840       -20      928           1007
##  7  2023     1    19     1932           1952       -20     2307           2301
##  8  2023     1    25      715            735       -20      916            950
##  9  2023     1    28     1435           1455       -20     1734           1805
## 10  2023     1     9      816            835       -19     1008           1040
## # ℹ 435,342 more rows
## # ℹ 11 more variables: arr_delay <dbl>, carrier <chr>, flight <int>,
## #   tailnum <chr>, origin <chr>, dest <chr>, air_time <dbl>, distance <dbl>,
## #   hour <dbl>, minute <dbl>, time_hour <dttm>

mutate

Pour créer de nouvelles colonnes à partir de colonnes existantes

head(select(flights, air_time), 5)

## # A tibble: 5 × 1
##   air_time
##      <dbl>
## 1      367
## 2      108
## 3      190
## 4      108
## 5       80

Création d’une colonne air_time_hours

flights <- mutate(flights, air_time_hours = air_time / 60)
head(select(flights, air_time, air_time_hours), 5)

## # A tibble: 5 × 2
##   air_time air_time_hours
##      <dbl>          <dbl>
## 1      367           6.12
## 2      108           1.8 
## 3      190           3.17
## 4      108           1.8 
## 5       80           1.33

Possible de créer plusieurs nouvelles colonnes d’un coup

flights <- mutate(
    flights,
    distance_km = distance / 0.62137,
    vitesse = distance_km / air_time_hours
)
select(flights, air_time_hours, distance_km, vitesse)

## # A tibble: 435,352 × 3
##    air_time_hours distance_km vitesse
##             <dbl>       <dbl>   <dbl>
##  1           6.12       4023.    658.
##  2           1.8        1223.    680.
##  3           3.17       2536.    801.
##  4           1.8        1024.    569.
##  5           1.33        785.    589.
##  6           2.57       1746.    680.
##  7           3.2        2536.    793.
##  8           1.98       1157.    583.
##  9           4.3        2253.    524.
## 10           2.62       1714.    655.
## # ℹ 435,342 more rows

Le pipe

Les actions vues précédemment peuvent être enchaînées grâce au pipe c’est à dire l’envoi du résultat de l’action précédente comme premier argument de l’action suivante. C’est l’opérateur %>% qui joue ce rôle

flights %>%
  filter(dep_delay >= 10 & dep_delay <= 15) %>%
  arrange(desc(dep_delay)) %>%
  select(flight, origin, dest, dep_delay, air_time_hours, distance_km, vitesse) %>%
  head(5)

## # A tibble: 5 × 7
##   flight origin dest  dep_delay air_time_hours distance_km vitesse
##    <int> <chr>  <chr>     <dbl>          <dbl>       <dbl>   <dbl>
## 1    786 EWR    RSW          15           2.73       1719.    629.
## 2    910 LGA    MSY          15           2.92       1904.    653.
## 3    908 EWR    MSY          15           2.88       1878.    651.
## 4    503 EWR    SFO          15           5.68       4128.    726.
## 5     27 JFK    SLC          15           4.97       3203.    645.

Attention ici les actions s’enchaînent (pipeline) mais rien n’est conservé. Si on veut placer le résultat dans un objet il faut commencer par objet <-

delays <- flights %>%
  filter(dep_delay >= 10 & dep_delay <= 15) %>%
  arrange(desc(dep_delay)) %>%
  select(flight, origin, dest, dep_delay, air_time_hours, distance_km, vitesse)

Cette fois le résultat du pipeline sera stocké dans l’objet delays

Pratique

• afficher les lignes 100 à 105 du dataset flights
• afficher les vols du mois de juillet 2013
• afficher les vols qui ont eu un retard (variable arr_delay) compris entre 5 et 15. Le résultat doit être affiché par ordre décroissant
• selectionner les colonnes name, lat et lon du dataset airports et renommer lat et lon en latitude et longitude
• dans le dataset airports, créer une nouvelle variable alt_metres qui contient l’altitude en mètres (conversion de pieds en mètres en divisant la variable alt par 3.28)

Grouper les données avec group_by

• permet de grouper les données puis d’appliquer une fonction à ces données groupées
• exemple: grouper les données par mois et affiche le premier index (combinaison group_by et slice)

flights %>% group_by(month) %>% slice(1)

## # A tibble: 12 × 22
## # Groups:   month [12]
##     year month   day dep_time sched_dep_time dep_delay arr_time sched_arr_time
##    <int> <int> <int>    <int>          <int>     <dbl>    <int>          <int>
##  1  2023     1     1        1           2038       203      328              3
##  2  2023     2     1      453            501        -8      934            943
##  3  2023     3     1        6           2259        67      104              5
##  4  2023     4     1        1           2205       116      111           2317
##  5  2023     5     1       10           2225       105      357            227
##  6  2023     6     1       39           1659       460      325           2000
##  7  2023     7     1        2           2256        66      401            300
##  8  2023     8     1        3           2035       208      224           2330
##  9  2023     9     1       12           2236        96      120           2359
## 10  2023    10     1       58           2250       128      302            118
## 11  2023    11     1        2           2340        22       57             54
## 12  2023    12     1      506            509        -3      807            823
## # ℹ 14 more variables: arr_delay <dbl>, carrier <chr>, flight <int>,
## #   tailnum <chr>, origin <chr>, dest <chr>, air_time <dbl>, distance <dbl>,
## #   hour <dbl>, minute <dbl>, time_hour <dttm>, air_time_hours <dbl>,
## #   distance_km <dbl>, vitesse <dbl>

Grouper les données par compagnie, calculer le delay moyen et le mettre dans une varaible mean_delay (combinaison group_by et mutate)

flights %>%
  group_by(carrier) %>%
  mutate(mean_delay = mean(dep_delay, na.rm=T)) %>%
  select(carrier, origin, dest, mean_delay)

## # A tibble: 435,352 × 4
## # Groups:   carrier [14]
##    carrier origin dest  mean_delay
##    <chr>   <chr>  <chr>      <dbl>
##  1 UA      EWR    SMF         17.6
##  2 DL      JFK    ATL         15.1
##  3 B6      JFK    BQN         23.8
##  4 B6      JFK    CHS         23.8
##  5 UA      EWR    DTW         17.6
##  6 AA      EWR    MIA         14.2
##  7 B6      JFK    BQN         23.8
##  8 AA      EWR    ORD         14.2
##  9 UA      EWR    IAH         17.6
## 10 NK      EWR    FLL         18.2
## # ℹ 435,342 more rows

Grouper les données par mois, et afficher les vols qui ont le délai maximum (combinaison group_by et filter)

flights %>%
  group_by(month) %>%
  filter(dep_delay == max(dep_delay, na.rm=T)) %>%
  select(carrier, month, origin, dest, dep_delay) %>%
  arrange(month)

## # A tibble: 12 × 5
## # Groups:   month [12]
##    carrier month origin dest  dep_delay
##    <chr>   <int> <chr>  <chr>     <dbl>
##  1 AA          1 JFK    STT        1201
##  2 AA          2 LGA    DFW        1665
##  3 AA          3 EWR    DFW        1557
##  4 AA          4 LGA    DCA        1742
##  5 UA          5 LGA    ORD        1312
##  6 AA          6 LGA    DFW        1353
##  7 AA          7 LGA    DFW        1413
##  8 F9          8 LGA    ATL        1161
##  9 AA          9 JFK    LAX        1369
## 10 AA         10 EWR    ORD        1781
## 11 DL         11 JFK    ATL        1074
## 12 AA         12 EWR    CLT        1813

Appliquer une fonction aux données: summarise

flights %>%
  summarise(
    retard_depart=mean(dep_delay, na.rm=T),
    retard_arrivee=mean(arr_delay, na.rm=T)
  )

## # A tibble: 1 × 2
##   retard_depart retard_arrivee
##           <dbl>          <dbl>
## 1          13.8           4.34

Utilisée avec group_by, summarise permet d’agréger une information sur des données groupées

On affiche le retard moyen de chaque compagnie et ordonne le résultat du plus grand retard moyen au plus petit

flights %>%
  group_by(carrier) %>%
  summarise(retard_depart=mean(dep_delay, na.rm=T)) %>%
  arrange(desc(retard_depart))

## # A tibble: 14 × 2
##    carrier retard_depart
##    <chr>           <dbl>
##  1 F9              35.7 
##  2 B6              23.8 
##  3 HA              22.9 
##  4 OO              19.8 
##  5 NK              18.2 
##  6 UA              17.6 
##  7 WN              16.1 
##  8 DL              15.1 
##  9 AA              14.2 
## 10 AS              12.0 
## 11 MQ              10.5 
## 12 9E               7.44
## 13 YX               4.21
## 14 G4               3.98

Si on veut connaître le nombre de lignes concernées, on passe l’argument n=n() dans l’appel de la fonction summarise

flights %>%
  group_by(carrier) %>%
  summarise(
    retard_depart=mean(dep_delay, na.rm=T),
    n=n()
  )%>%
  arrange(desc(retard_depart))

## # A tibble: 14 × 3
##    carrier retard_depart     n
##    <chr>           <dbl> <int>
##  1 F9              35.7   1286
##  2 B6              23.8  66169
##  3 HA              22.9    366
##  4 OO              19.8   6432
##  5 NK              18.2  15189
##  6 UA              17.6  79641
##  7 WN              16.1  12385
##  8 DL              15.1  61562
##  9 AA              14.2  40525
## 10 AS              12.0   7843
## 11 MQ              10.5    357
## 12 9E               7.44 54141
## 13 YX               4.21 88785
## 14 G4               3.98   671

Possible de grouper les données selon plusieurs variables
On groupe les données par mois et par destination et on compte le nombre de vols

flights %>% 
  group_by(month, dest) %>%
  summarise(n=n()) %>%
  arrange(desc(n))

## `summarise()` has grouped output by 'month'. You can override using the
## `.groups` argument.

## # A tibble: 1,274 × 3
## # Groups:   month [12]
##    month dest      n
##    <int> <chr> <int>
##  1     3 BOS    2008
##  2     5 BOS    1798
##  3     4 BOS    1768
##  4     2 BOS    1745
##  5     1 BOS    1693
##  6     3 ORD    1654
##  7     3 MCO    1624
##  8     1 ORD    1603
##  9     5 ORD    1603
## 10    12 MCO    1585
## # ℹ 1,264 more rows

On peut utiliser les variables générées par summarise pour créer une nouvelle variable

flights %>% 
  group_by(month, dest) %>%
  summarise(n=n()) %>%
  mutate(pourcentage=n/sum(n) * 100) %>%
  arrange(desc(pourcentage))

## `summarise()` has grouped output by 'month'. You can override using the
## `.groups` argument.

## # A tibble: 1,274 × 4
## # Groups:   month [12]
##    month dest      n pourcentage
##    <int> <chr> <int>       <dbl>
##  1     3 BOS    2008        5.08
##  2     2 BOS    1745        5.02
##  3    12 MCO    1585        4.75
##  4     4 BOS    1768        4.72
##  5     1 BOS    1693        4.70
##  6     5 BOS    1798        4.64
##  7     1 ORD    1603        4.45
##  8    11 ORD    1514        4.39
##  9     2 ORD    1519        4.37
## 10    11 MCO    1496        4.33
## # ℹ 1,264 more rows

Pratique

• afficher le nombre de vols par mois et afficher le résultat par ordre croissant du nombre de vols
• calculer le nombre de vols à destination de Los Angeles (code LAX) pour chaque mois de l’année
• calculer le nombre de vols par mois et par destination et afficher uniquement le max
• calculer le nombre de vols par mois et créer une colonne qui affiche le pourcentage que cela représente sur les vols annuels

Graphiques avec ggplot2

dataset pour travailler

rp <- rp2012 %>% filter(departement %in% c("Oise", "Rhône", "Hauts-de-Seine", "Lozère", "Bouches-du-Rhône"))
head(rp)

## # A tibble: 6 × 60
##   code_insee commune     code_region region code_departement departement pop_tot
##   <chr>      <chr>       <chr>       <chr>  <chr>            <chr>         <int>
## 1 13001      Aix-en-Pro… 93          Prove… 13               Bouches-du…  141148
## 2 13002      Allauch     93          Prove… 13               Bouches-du…   20690
## 3 13003      Alleins     93          Prove… 13               Bouches-du…    2428
## 4 13004      Arles       93          Prove… 13               Bouches-du…   52439
## 5 13005      Aubagne     93          Prove… 13               Bouches-du…   45243
## 6 13007      Auriol      93          Prove… 13               Bouches-du…   11621
## # ℹ 53 more variables: pop_cl <fct>, pop_0_14 <int>, pop_15_29 <int>,
## #   pop_18_24 <int>, pop_75p <int>, pop_femmes <int>, pop_act_15p <int>,
## #   pop_chom <int>, pop_agric <int>, pop_indep <int>, pop_cadres <int>,
## #   pop_interm <int>, pop_empl <int>, pop_ouvr <int>, pop_scol_18_24 <int>,
## #   pop_non_scol_15p <int>, pop_dipl_aucun <int>, pop_dipl_bepc <int>,
## #   pop_dipl_capbep <int>, pop_dipl_bac <int>, pop_dipl_bac2 <int>,
## #   pop_dipl_sup <int>, log_rp <int>, log_proprio <int>, log_loc <int>, …

Histogramme

• d’abord initialiser la fonction graphique en lui passant le dataset utilisé: ggplot(data=rp)
• ensuite on ajoute les graphiques (geom) avec un “+”, ici geom_histogram
• dans le geom on ajoute les arguments demandés, notamment dans la fonction aes()

ggplot(data=rp) +
  geom_histogram(aes(x=cadres))

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Nuage de points

il faut passer 2 arguments à la fonction aes (aesthetic), la variable x et la variable y

ggplot(data=rp) +
  geom_point(aes(x=dipl_sup, y=cadres))

Pour personnaliser le graphique on peut passer d’autres arguments à geom_point en plus de la fonction aes

ggplot(data=rp) +
  geom_point(aes(x=dipl_sup, y=cadres), col="blue", size=2)

Boxplot

ggplot(data=rp) +
  geom_boxplot(aes(x=departement, y=maison))

L’argument varwidth permet que la largeur de la box soit proportionnelle au nombre d’observations

ggplot(data=rp) +
  geom_boxplot(aes(x=departement, y=maison), varwidth=T)

Graphique en barres

ggplot(data=rp) +
  geom_bar(aes(x=departement))

Si on veut ajouter un “y” que l’on a calculé soit même alors il faut ajouter l’argument stat=“identity” dans l’appel de la fonction geom_bar
Exemple: on calcule le pourcentage que représente le nombre de communes de chaque département par rapport au total, en utilisant les connaissances acquises

rp %>% group_by(departement) %>% summarise(n=n()) %>% mutate(pourcent=n/sum(n) * 100) %>%
  ggplot() + geom_bar(aes(x=departement, y=pourcent), stat="identity")

Courbes d’évolution dans le temps

chargement d’un jeu de données du package ggplot2 : economics

data("economics")
economics

chom <- economics %>% select(date, unemploy)
ggplot(chom) + 
  geom_line(aes(x=date, y=unemploy))

si on souhaite plusieurs courbes on fait plusieurs appels à geom_line (avec “+” entre chaque)

Principe du mappage

• Associer un des attributs du graphique (couleur, taille) à une variable
• Permet de rendre visuelle une troisième dimension (voire une quatrième)
• Le mappage se fait dans l’appel de la fonction aes()

ggplot(rp) +
  geom_point(aes(x = dipl_sup, y = cadres, color = departement))

on ajoute une 4ème dimension à travers la taille des points, qui représente la population totale de la commune

ggplot(rp) +
  geom_point(aes(x=dipl_sup, y=cadres, color=departement, size=pop_tot))

Possible de paramétrer un mappage avec scale_

ggplot(rp) +
  geom_point(aes(x=dipl_sup, y=cadres, color=departement, size=pop_tot)) +
  scale_size("Population", breaks=c(0, 5000, 10000, 50000, 100000), range=c(0, 10))

Ajouter une droite de régression dans un nuage de points avec geom_smooth

Si on ajoute plusieurs geom qui utilisent les mêmes arguments dans la fonction aes, on peut passer la fonction aes() dans l’appel de la fonction ggplot

ggplot(rp, aes(x=dipl_sup, y=cadres)) +
  geom_point() + geom_smooth(method ="lm", color="red")

## `geom_smooth()` using formula = 'y ~ x'

Principe du faceting

• répéter un graphique plusieurs fois selon une ou plusieurs variables

ggplot(rp) +
  geom_histogram(aes(x=cadres)) +
  facet_wrap(vars(departement))

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Personnalisation des axes

scale_x_continuous et scale_y_continuous pour les axes x et y de données quantitatives

ggplot(data=rp) +
  geom_point(aes(x=dipl_sup, y=cadres)) +
  scale_x_continuous("Diplômes supérieurs", limits=c(0, 60)) +
  scale_y_continuous("Cadres", limits=c(0, 60))

scale_x_discrete et scale_y_discrete pour les variables qualitatives

ggplot(rp) +
  geom_bar(aes(x = departement)) +
  scale_x_discrete("Département")

Rmq: il existe un package, esquisse, qui propose ensuite un Addins pour faire les graphiques avec une interface graphique: install.packages(“esquisse”), puis library(esquisse) et ensuite Addins/ggplot2 builder

Pratique

Avec le dataset rp2012

• un nuage de points x=dipl_aucun et y=ouvr, avec titre des axes, et points en bleu
• même graphique que précédemment mais avec un mappage de couleur selon la région
• boxplot: distribution de la variable proprio par rapport à la taille de departement (pop_cl)
• nuage de points avec en x=dipl_sup et y=cadres avec un graphique par région (facet_wrap)
• boxplot: distribution du pourcentage moyen de propriétaires par département au sein de chaque région

• un ouvrage complet sur l’analyse de données avec R ici
• le cours de Julien Barnier accessible en ligne ici
• le site https://www.tidyverse.org/
• le site R for data science
• le site de ggplot2 http://ggplot2.tidyverse.org/reference/
• des exemples de graphiques avec ggplot2 : http://www.cookbook-r.com/Graphs/